This PR contains the main algorithm for the new decimal <--> floating conversion code. This algorithm was written to address the precision issues described here.

Summary

• The new algorithm is more accurate than the previous code, but it is also far more complex.
• It can perform conversions that were not even possible in the old code due to overflow (decimal32/64/128 conversions only worked for scale factors up to 10^9/18/38, respectively). Now the entire floating-point range is convertible, including denormals.
• This new algorithm is significantly faster in some parts of the conversion phase-space, and in some parts slightly slower.

Previous PR's

These contain the supporting parts of this work:

• Explicit conversion PR
• Benchmarking PR
• Powers-of-10 PR
• Utilities PR. These utilities are updated here to support denormals.

Algorithm Outline

We convert floating -> (integer) decimal by:

• Extract the floating-point mantissa (converted to integer) and power-of-2
• For float we use a uint64 to contain our data during the below shifting/scaling, for double uint128_t
• In this shifting integer, we alternately apply the extracted powers-of-2 (bit-shifts, until they're all used) and scale-factor powers-of-10 (multiply/divide) as needed to reach the desired scale factor.

Decimal -> floating is just the reverse operation.

Supplemental Changes

• Testing: Add decimal128, add precise-conversion tests. Remove kludges due to inaccurate conversions. Add test for zeroes.
• Benchmarking: Enable regions of conversion phase-space for benchmarking that were not possible in the old algorithm.
• Unary: Cleanup by using CUDF_ENABLE_IF. Call new conversion code for base-10 fixed-point.

Performance for various conversions/input-ranges

• Note: F32/F64 is float/double

New algorithm is FASTER by:

• F64 --> decimal64: 60% for E8 --> E15
• F64 --> decimal128: 13% for E-8 --> E-15
• F64 --> decimal128: 22% for E8 --> E15
• F64 --> decimal128: 27% for E31 --> E38
• decimal32 --> F64: 18% for E-3 --> E4
• decimal64 --> F64: 27% for E-14 --> E-7
• decimal64 --> F64: 17% for E-3 --> E4
• decimal128 --> F64: 21% for E-14 --> E-7
• decimal128 --> F64: 11% for E-3 --> E4
• decimal128 --> F64: 13% for E31 --> E38

New algorithm is SLOWER by:

• F32 --> decimal32: 3% for E-3 --> E4
• F32 --> decimal64: 2% for E-14 --> E14
• F64 --> decimal32: 3% for E-3 --> E4
• decimal32 --> F32: 5% for E-3 --> E4
• decimal128 --> F64: 36% for E-37 --> E-30

Other kernels:

• The PYMOD binary-op benchmark is 7% slower.

Performance discussion

• Many conversions have identical speed, indicating these algorithms are often fast and we are instead bottlenecked on overheads such as getting the input to the gpu in the first place.
• F64 conversions are often much faster than the old algorithm as the new algorithm completely avoids the FP64 pipeline. Other than the cast to double itself, all of the operations are on integers. Thus we don't have threads competing with each other and taking turns for access to the floating-point cores.
• The conversions are slightly slower for floats with powers-of-10 near zero. Presumably this is due to code overhead for e.g., handling a large range of inputs, UB-checks for bit shifts, branches for denormals, etc.
• The conversion is slower for decimal128 conversions with very small exponents, which requires several large divisions (128bit divided by 64bit).
• The PYMOD kernel is slower due to register pressure from the introduction of the new division routines in the earlier PR. Even though this benchmark does not perform decimal <--> floating conversions, it gets hit because of inlined template code in the kernel increasing the code/register pressure.

Checklist

• I am familiar with the Contributing Guidelines.
• New or existing tests cover these changes.
• The documentation is up to date with these changes.